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Abstract 


We report the discovery of a rare new form of long-term radio variability in the light curves of active galaxies (AG) 
—symmetric achromatic variability (SAV)—a pair of opposed and strongly skewed peaks in the radio flux density 
observed over a broad frequency range. We propose that SAV arises through gravitational milli-lensing when 
relativistically moving features in AG jets move through gravitational lensing caustics created by 102-106 Mo 
subhalo condensates or black holes located within intervening galaxies. The lower end of this mass range has been 
inaccessible with previous gravitational lensing techniques. This new interpretation of some AG variability can 
easily be tested and if it passes these tests, will enable a new and powerful probe of cosmological matter 
distribution on these intermediate-mass scales, as well as provide, for the first time, micro-arcsecond resolution of 
the nuclei of AG—a factor of 30-100 greater resolution than is possible with ground-based millimeter very-long- 


baseline interferometry. 


Key words: BL Lacertae objects: individual (J1415+ 1320) — gravitational lensing: strong — radio continuum: 
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1. Introduction 


We report a hitherto unrecognized form of radio variability 
in active galaxies (AG)—symmetric achromatic variability 
(SAV)—which is time-symmetric and achromatic from 15 to 
234 GHz. SAV was first noticed in an unusual year-long 
symmetric U-shaped feature that appeared in the 15 GHz light 
curve of the BL Lac object J14154-1320 (PKS 1413--135), in 
2009-2010 (see Figure 1). The radio source lies either within or 
behind a spiral galaxy at redshift z — 0.247 (Beichman et al. 
1981; McHardy et al. 1991; Carilli et al. 1992). The U-shaped 
symmetry is similar to that observed in extreme scattering 
events (ESEs; Fiedler et al. 1987; Pushkarev et al. 2013) caused 
by interstellar plasma structures, but the achromaticity, among 
other factors, can be shown to rule out ESEs as the origin of 
SAV (Vedantham et al. 2017). A pair of forward-reverse 
shocks within the source may create achromatic U-shaped 
events (Kobayashi & Zhang 2003), but these are not expected 
to yield the high degree of symmetry observed in J14154-1320. 

We propose that SAV is caused by the modulation of the 
flux density of a luminal or superluminal (Begelman et al. 
1984) compact emission region as it traverses the magnification 
pattern cast by an intervening gravitational lens. The proposed 
lensing mechanism is akin to microlensing of stars by stellar- 
mass lenses (Alcock et al. 1993), but the lenses would need to 


10 We use “SAV” to denote symmetric achromatic variability, symmetric 
achromatic variations, and symmetric achromatic variable, relying on the 
context to make the usage clear. 


be in the milli-lensing mass range ~10°-10°M,, with a 
projected surface mass density of >10*M. pe *. Potential lens 
candidates with these properties include the dense cores of 
globular clusters and molecular clouds, and massive black 
holes. The lenses likely consist of multiple components that are 
projected close to the line of sight to the source, but are not 
necessarily gravitationally bound to each other. 

The intermediate-mass range has not been accessible by 
traditional centimeter-wavelength very-long-baseline interfero- 
metry (VLBI) searches for multiple images (Wilkinson et al. 
2001) or time delays (Fassnacht et al. 1999), and although the 
upper half of this mass range is now being probed by 
millimeter wavelength VLBI (Marscher & Jorstad 2010; 
Brinkerink et al. 2016), the lower half is still barely within 
reach. SAV may thus provide a powerful new method of 
exploiting gravitational lensing in addition to the two 
traditional methods employing multiple images and time 
delays. We note that our gravitational milli-lensing hypothesis 
leads to a number of testable predictions and is correspondingly 
fragile. We regard this as a strength of the hypothesis—it 
should be easy to disprove. 

The paper is organized as follows. Section 2 describes the 
Owens Valley Radio Observatory (OVRO) blazar monitoring 
survey in which SAV was first identified, and presents 
observational details of multifrequency radio light curves of 
J14154-1320. Section 3 describes the SAV observed in J1415 
4-1320 and the rate of U-shaped events in the OVRO survey. In 
Section 4 we summarize the known radio and optical properties 
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Figure 1. Multifrequency light curves of J1415--1320 showing symmetric achromatic variability (SAV). The roughly 1-year-long U-shaped SAV events in 2009 and 


2014 are highlighted. 


of J14154-1320. In Section 5 we present simple lens models 
that can account for the observed SAV features in J1415-4- 1320 
as an "existence proof." Though the light curves cannot be 
inverted to obtain a unique lens model, we exploit some generic 
properties of gravitational lenses to place constraints on the 
lens mass and projected surface density. Finally, in Section 6, 
we discuss the general and J1415+1320-specific implications 
of the lensing hypothesis, and make falsifiable predictions to 
test the lensing hypothesis. 


2. The Observations 


We report here on four series of radio flux density- 
monitoring observations, carried out at 15 GHz on the 40m 
telescope at OVRO, at 22 and 37 GHz on the 13.7 m telescope 
of the Metsáhovi Radio Observatory (MRO), at 100 GHz on 
the 6 x 10.4 m OVRO millimeter array, and at 234 GHz on the 
8 x 6m Sub-Millimeter Array (SMA) of the Smithsonian 
Astrophysical Observatory (SAO). 


2.1. OVRO 15 GHz 
2.1.1. The Sample 


Since 2008 the OVRO 40 m telescope blazar monitoring 
survey has accumulated ~12,000 object-years of observations 
and is by far the largest and most sensitive radio monitoring 
survey of blazars undertaken thus far. The survey is being 
carried out in support of the Fermi Gamma-Ray Space 
Telescope and the sample comprises mostly candidate 
gamma-ray blazars. This program started as a survey of 
15 GHz radio variability of a complete sample of 1158 AG 
(Richards et al. 2011), which has been augmented by the 
addition of several hundred Fermi-detected blazars that are not 
in the original sample, samples of galactic micro-quasars, 
NLS1 galaxies, and the sample of sources detected at TeV 
energies,!! so that our whole sample now contains ~1830 
objects each observed about twice per week, hardware 
problems and weather permitting. 


11 http: //tevcat.uchicago.edu/ 


The complete sample comprising 1158 objects was selected 
from a uniform all-sky survey of bright quasars in the 
Candidate Gamma-Ray Blazar Survey (CGRaBS; Healey 
et al. 2008). The CGRaBS sample was selected from objects 
in the Cosmic Lens All Sky Survey (Myers et al. 2003), 
together with observations at both 4.8 and 8.4GHz. Over 
11,000 objects brighter than 65 mJy at 4.8 GHz, and with 
spectral index a > —0.5 (S x rv?) were selected initially. Since 
the 40 m program is in support of Fermi-GST, we restricted the 
sample to EGRET-like blazars (Hartman et al. 1999; Mattox 
et al. 2001), and thus applied a figure of merit based on 3EG 
blazar objects in CLASS, which yielded 1625 objects, of which 
1158 have declinations >—20° and can be monitored by the 
40m telescope. We concentrate here on this statistically 
complete sample, in which we find that 981 objects have a 
high enough signal-to-noise ratio to identify symmetrical 
U-shaped features in their light curves. 


2.1.2. The Observations 


The 15 GHz observations were made with the 40 m telescope 
of the OVRO in California (3772314 N, 11872827 W). The 
40 m telescope uses off-axis dual-beam optics and a cryogenic 
receiver with a 15 GHz center frequency and 3 GHz bandwidth. 
The two sky beams are Dicke-switched using the off-source 
beam as a reference. The source is alternated between the two 
beams in an ON-ON fashion to remove atmospheric and 
ground contamination. The on-source observing time is 32s 
and the 3c detection limit is 12 mJy in good weather. 

In 2014 May, a new pseudo-correlation receiver was 
installed on the 40 m telescope, and the gain variations have 
since been corrected using a 180? phase switch instead of a 
Dicke switch. The performance of the new receiver is very 
similar to the old one and no discontinuity is seen in the light 
curves. Calibration is achieved using a temperature-stable 
diode noise source to remove receiver gain drifts, and the flux 
density scale is derived from observations of 3C 286 assuming 
a value of 3.44 Jy at 15 GHz (Baars et al. 1977). The systematic 
uncertainty of about 596 in the flux density scale is not included 
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Figure 2. Same as Figure 1 but with the entire Metsähovi monitoring data at 22 and 37 GHz since 1989, and the OVRO MMA data at 100 GHz. Two more achromatic 


U-shaped events are seen in 1993 and 2000. 


in the error bars. Further details of the reduction and calibration 
procedures may be found in Richards et al. (2011). 


2.2. Metsähovi 22 and 37 GHz 


The 22 and 37 GHz observations were made with the 13.7 m 
diameter Aalto University Metsáhovi radio telescope, which is 
a radome-enclosed Cassegrain antenna in Finland (607218 N, 
247393 E). The measurements were made with 1 GHz-band 
dual-beam receivers centered at 22.2 and 36.8 GHz, respec- 
tively. The high electron mobility transistor front end operates 
at room temperature. The observations are Dicke-switched 
ON-ON observations, alternating the source and an adjacent 
patch of sky in each feed horn. The typical integration time for 
one flux-density measurement is 1200-1600 s. The detection 
limit is about 0.2 Jy under optimal conditions. Data points with 
a signal-to-noise ratio «4 are considered non-detections. The 
flux-density scale is set by observations of the HII region 
DR21. Sources NGC 7027, 3C274, and 3C 84 are used as 
secondary calibrators. A detailed description of the data 
reduction and analysis is given in Teraesranta et al. (1998). 
The error estimate in the flux density includes the contribution 
from the measurement rms and the uncertainty of the absolute 
calibration. 


2.3. OVRO 100 GHz 


Flux densities in the 3 mm band (100 GHz) were obtained in 
the period 1992-2003 from the flux density history catalog of 
the OVRO Millimeter Array (MMA), a six-element millimeter- 
wave radio telescope array located at the main OVRO site until 
2005. The observations were primarily made during dedicated 
calibration sessions to track the flux density of millimeter-loud 
blazars and other point sources. The flux-density scale was 
calibrated using contemporaneous measurements of Uranus 
and/or Neptune. 


2.4. SMA 234 GHz 


Flux densities were measured at the Submillimeter Array 
(SMA), an 8-element interferometer located near the summit of 


Maunakea (Hawaii). The measurements were made at frequen- 
cies ranging from 221 to 297 GHz; we use the mean frequency 
of 234 GHz in our analysis. J1415+1320 is included in an 
ongoing monitoring program at the SMA to determine the flux 
densities of compact extragalactic radio sources that can be 
used as calibrators at millimeter and submillimeter wavelengths 
(Gurwell et al. 2007). Potential calibrators are from time to time 
observed for 3-5 minutes, and calibrated against known 
standards, typically solar system objects (Titan, Uranus, 
Neptune, or Callisto). Data from this program are updated 
regularly and are available at the SMA website.” 


3. Symmetric Achromatic Variability (SAV) 


In 2009, we observed a roughly 1-year-long highly 
symmetric U-shaped feature in the OVRO light curve of 
J14154-1320 and a similar feature was seen in 2014. Figure 1 
shows the 15 GHz (OVRO) 37GHz (Metsáhovi) and 
234 GHz (SAO) light curves. The U-shaped features are highly 
achromatic from 15 to 234 GHz, although there are small 
deviations from perfect achromaticism. We therefore refer to 
this type of variability as “SAV.” The object has also been 
monitored at 22 and 37 GHz at Metsáhovi since 1989, and two 
more U-shaped features were observed in 1993 and 2000, 
which were also observed on the OVRO MMA at 100 GHz 
(Figure 2). Thus, these two additional U-shaped features are 
also achromatic over the whole observed range—here from 22 
to 100 GHz. 

The details of the profile of the U-shaped features are seen 
most clearly in the high-cadence, high-sensitivity OVRO 
observations from 2008 to 2016 (Figure 1). These show a 
distinctive mirror symmetry composed of a pair of spikes with 
time-reversed profiles. As described in Section 5 in the context 
of gravitational lensing, we can distinguish two types of 
U-shaped events: a "volcano" type in which the first spike 
displays a slow rise and a fast decline and the second spike 
shows a fast rise and slow decline (SRFD-FRSD); and a 
"crater" type in which a fast rise, slow decline is followed by a 


12 http: / [smal.sma.hawaii.edu/callist/callist.html 


THE ASTROPHYSICAL JOURNAL, 845:89 (16pp), 2017 August 20 


slow rise, fast decline (FRSD-SRFD). All four events seen in 
J1415+1320 are of the volcano-type. We see from Figures 1 
and 2 that an eight-year period of low intrinsic activity has 
enabled us to detect SAV clearly in J1415+1320 between 2008 
and 2016. 


3.1. Incidence of U-shaped Events 


Many sources in the complete sample often show “flares” 
that can have either SRFD or FRSD shape (Hovatta et al. 
2008). To test if the observed U-shaped events are due to a 
chance juxtaposition of two flares, we manually examined our 
complete sample to derive an approximate estimate for the 
incidence of U-shaped events of both crater and volcano-type. 
981 objects in our sample had light curves of sufficient quality 
to identify U-shaped events. Among these, we identified 10 
U-shaped events in 7 objects (see Figure 9). The probability of 
a U-shaped feature in 1 of our 981 8 year light curves is 
therefore ~1%. Even among this small sample, J1415+1320 
stands out in view of the clarity and isolation of the two 
U-shaped features relative to the rest of the light curves. 

If such U-shaped features are distributed randomly among all 
the objects, the probability of three or more occurrences in one 
object in 27 years (J1415+1320) is ~4 x 10-3. Though our 
manual selection criterion is admittedly crude, the extremely 
low probability suggests that the repeating U-shaped features in 
J1415+1320 are not random events but are an unusual intrinsic 
feature of the source itself, due to propagation effects along the 
lines of sight, or that are a combination of both. Besides this 
statistical argument, J14154-320 has other unusual properties 
that have previously led to suggestions that it is being lensed. 
We explore them in Section 4 below. 


4. Radio and Optical Properties of J1415+1320 


J1415+1320 has a two-sided jet and exhibits up to eight 
mas-scale components extending over ~60 mas, some of which 
move outward with speeds ranging from 0.06 to 0.3 mas yr! 
and apparent transverse speeds up to 1.5c (Kellermann et al. 
1998; Perlman et al. 2002; Gugliucci et al. 2005; Lister et al. 
2009, 2016). Based on its radio morphology, J1415+1320 is 
classified as a compact symmetric object (CSO), with radio 
lobes straddling a bright core (Conway et al. 1992; Wilkinson 
et al. 1994; Readhead et al. 1996). Although the cores in CSOs 
are bright (Conway et al. 1992) and variable, the core of J1415 
+1320 is exceptional—with a core fraction 70.5 (Lister et al. 
2009) it is both far brighter and far more variable than most 
CSO cores in radio galaxies, which have core fractions in the 
range ^0.004—0.06 (Taylor et al. 1996). In addition, J1415 
+1320 is identified with a spiral galaxy at z — 0.247 viewed 
edge-on (Beichman et al. 1981; McHardy et al. 1991; Carilli 
et al. 1992), whereas BL-Lac objects are rarely, if ever, hosted 
by spiral galaxies. The radio core is offset from the optical 
center of the spiral galaxy by 0.013 arcsec (50 pc) (Perlman 
et al. 2002). These unusual properties, among others, have led 
to suggestions that the radio source is an unrelated background 
object that is being lensed by matter in the spiral galaxy (Carilli 
et al. 1992; Stocke et al. 1992). We note that based on the 
number densities of radio sources and intervening galaxies, the 
probability of such a chance alignment is small (20.006; 
details of calculation in Appendix D). However, this prob- 
ability is based on a posteriori statistics and as such the 
background source hypothesis cannot be ruled out. If the radio 
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source is a background object, then due to the lack of multiple 
images on arcsecond scales (due to lensing by the bulge of the 
spiral), the redshift of the radio source is bounded by z € 0.5 
(Lamer et al. 1999). 


4.1. Optical Spectrum of J1415+1320 


We obtained an optical spectrum of the spiral galaxy with 
the low-resolution imaging spectrometer (LRIS; Oke et al. 
1995) on the Keck telescope. The data were obtained on 
MJD 57488.43 with a total exposure time of about 1 hour. The 
lpipe? software was used for bias subtraction, flat-fielding, 
cosmic-ray rejection, and sky-line subtraction. We then 
obtained a wavelength solution using arc-lamps and sky lines, 
optimally extracted the spectrum, and calibrated the flux- 
density and atmospheric opacity. The resulting spectrum is 
shown in Figure 3. We detected several emission and 
absorption lines with high significance, all of which can be 
attributed to the spiral galaxy at a redshift consistent with that 
obtained from 21cm absorption (Carilli et al. 1992). The 
spectrum also contains some lines (not annotated in Figure 3) 
detected at low significance that are consistent with z — 0 and 
likely caused by Galactic absorption. We find no clear evidence 
for an intervening galaxy. We also did not detect broad lines 
from the AGN, owing to severe extinction by foreground dust- 
lanes in the spiral galaxy. As such, the redshift of the AGN 
could not be independently determined. 


5. SAVs from Gravitational Lensing 


As shown in Section 3, the statistics of occurrence of 
symmetric U-shaped events makes them unlikely to be caused 
by random intrinsic fluctuations. In a companion paper 
(Vedantham et al. 2017), we have studied in detail the 
possibility that SAV could be caused by ESEs and shown that 
this explanation for the SAV phenomenon can be ruled out. In 
this section, we show that the SAVs we observe in J1415 
+1320 can be explained naturally and economically by 
gravitational lensing. 

Gravitational lensing is frequency-independent and accounts 
for the achromatic nature of SAVs. Relatively simple lens 
configurations yield symmetric magnification patterns. Time- 
dependent magnification occurs as a source moves over a 
magnification pattern cast by a foreground lens mass. The most 
prominent magnification features expected from gravitational 
lensing are fold caustics (Blandford & Narayan 1986, 1992; 
Blandford & Kovner 1988). Fold caustics are curves in the 
source plane where the point source magnification diverges." 
Rapid brightness changes ensue as a source crosses a fold 
caustic, and U-shaped features are a result of a closely spaced 
pair of fold crossings. Such U-shaped events, due to their 
symmetry, are a telltale sign of lensing, and a powerful 
discriminant between intrinsic and lensing-induced variability 
(Schneider et al. 1992). 

The generic magnification pattern in gravitational lensing 
comprises fold curves meeting at cusp points. When a source 
crosses into a closed fold contour at time t., the flux variation 
shows a rapid rise with a relatively slow decay that follows a 
|t — t.[-!/2 law (Schneider et al. 1992). This yields an FRSD 
feature. Source motion in the opposite direction naturally yields 


m http: //www.astro.caltech.edu / - dperley /programs/lpipe.html 


14 If we allow the distance of the source to vary, they should be seen as 
surfaces. 
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Figure 3. Continuum normalized LRIS spectrum of the putative host of J1415+1320. Prominent lines are marked for a redshift of z = 0.24675, which agrees with that 
determined from H I 21 cm absorption. There is no evidence for emission from an intervening redshift. The artifacts near 6877 and 7613 A are due to imperfect 


subtraction of sky lines. 


an SRFD curve. Crater symmetries (FRSD-SRFD) occur when 
a source crosses two folds that meet at a nearby cusp. Volcano 
symmetries (SRFD-FRSD) occur when a source crosses out of 
a closed fold contour and subsequently crosses back into it. 
Successive U-shaped SAV features may or may not come from 
the same source component; if they come from successive 
components following a similar track, then we expect the time 
series to be similar, although the transverse speed may be 
different, leading to a temporal scaling of the flux density 
variations. 


5.1. An Existence Proof 


To establish a “proof of existence" for the gravitational 
lensing hypothesis, we explored various commonly encoun- 
tered lens configurations that can yield the basic unit of SAV— 
a U-shaped event with volcano (SRFD-FRSD) symmetry—as 
a result of rectilinear motion of the source. Lenses with a high 
degree of circular symmetry typically yield a four-cusped 
diamond-shaped caustic configuration (called an astroid). 
Barring contrived source trajectories, such a lens yields 
crater-type light curves. The simplest asymmetric lens config- 
uration is a binary lens (the binary need not be gravitationally 
bound). Binary lenses typically yield a six-cusped caustic 
pattern that has sufficient complexity to generate a large range 
of light curve features including both crater and volcano- 
type SAVs. 

We simulated several binary lenses by varying the binary 
separation d and mass ratio q < 1 and found that a range of 
binary lens parameters can yield the observed U-shaped 
features for a variety of source trajectories (details of the 
simulations are presented in Appendix B). We show two 
examples in Figure 4. The upper two panels (left) show the 
magnification pattern of a binary lens with a mass ratio of 
q= 0.01 and separation of d = 0g (0g is the Einstein radius of 
the more massive object). Since such binary condensates 


typically exist in the presence of a smoothly varying 
background mass distribution (such as in a galaxy), we added 
a background convergence « = 0.1 and shear y = —0.05. These 
values are not very constraining, but are chosen to be "sub- 
critical, in that the external mass distribution cannot strongly 
lens the source independently. An approximately south-north 
source trajectory across the lens (black arrows) generates a 
single SRFD-FRSD U-shaped event. Hence two compact 
source components can generate the entire 8-year light curve to 
good accuracy. Note that we have removed a linear trend from 
the 15 GHz data points to account for any intrinsic variability 
of the lensed components or that of unlensed source 
components away from the caustics. 

The bottom panel of Figure 4 shows an alternative binary 
model that can reproduce both U-shaped features with only one 
moving radio component rather than two. However, this model 
requires the source to “graze” the caustic pattern, making the 
set of plausible source trajectories somewhat restrictive. 

The above fits were obtained by manually varying the lens 
and source parameters. In particular, the crossing timescale, 
and hence the source speed, are largely set by the relative 
timing of the four peaks that comprise the two U-shaped events 
in 2009 and 2014. The source components’ sizes and flux 
densities are largely determined by the peak-magnification at 
the caustics. These two factors are somewhat degenerate, 
because the magnification near a fold caustic is expected to 
follow u = ap™!/?, where p = o,/Üg where the source is 
modeled as a Gaussian with standard deviation o;, and a = 1 
depends on the precise nature of the Fermat potential at the 
fold. Hence a weaker smaller source can generate a caustic 
peak with a similar flux density as a stronger larger source. 
However, by exploring a large parameter range, we found that 
the source components must be smaller than ~0.03 Og to 
account for the precise light curve profile around the caustic 
crossing. For the family of lens models considered here 
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Figure 4. Fits to the OVRO 15 GHz data for two lens models. A linear function of time has been subtracted from the data. The left panels show the magnification 
pattern, and the right panels show the model light curves. The top two rows have the same lens model but different trajectories for the two source components (red and 
green). The right panels show the light curves for the individual components (green and red) and their sum (black). The lens model in the bottommost row requires 
only one source component but the range of plausible trajectories is limited. Lens and source parameter values are shown in the embedded text. 


we found 2.5 $a 354.25. Table 1 summarizes the parameter 
values for two lens models presented in Figure 4. 

A hallmark of gravitational lensing is its achromaticity, 
although we expect some chromaticity due to variation of 


source size and centroid position with frequency. Such 
variation is expected in the emission from relativistic jets due 
to a combination of the evolution of the emitting population of 
relativistic electrons with time and location in the jet, and 
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Table 1 
Summary of Lens and Source Parameters Used in the Model Fits from Figure 4 
q d/0, K Y tg (years) Og (uas) v/c p—23,/0g 
0.25 12 0.38 0.1 7.54 81 0.33 0.01 x 0.01 
0.01 1.0 0.1 —0.05 6.31 87.3 0.42 0.04 x 0.01 


Note. The columns represent the binary mass ratio q, binary separation d, external convergence x, external shear y, Einstein crossing time fg, Einstein radius Og, 


transverse speed of the source v, , and source size p. The Einstein angle assumes a source brightness temperature of T, = 10! K, and the transverse speed values are 
based on a nominal source redshift of z, — 0.5 (cosmological time dilation included). 


synchrotron opacity effects. To test if the binary lens models 
from Figure 4 reproduce the 37 and 234 GHz observations, we 
assumed that (i) the source component size varies with 
frequency v as v-!, (ii) the flux densities of the unlensed 
source components have a spectral index of o — —0.7, and 
(iii) the lensed components have a flatter spectral index 
of a=—0.1. With these assumptions, we predicted the 
corresponding light curves at 37 and 234 GHz where data are 
available. Figure 5 shows the multifrequency data and the 
expected light curves for the lens model shown in the upper 
row of Figure 4. Despite simplistic assumptions on the source 
properties, the fits are quite good over a factor of 15 in 
wavelength. 

The data show small but significant departures from both 
model light curves, which may be due to intrinsic variability or 
due to our simplistic assumptions of both the source and lens. 
For instance, there is evidence for additional symmetries in the 
light curves (see Appendix C), which may be accounted for by 
lenses with a larger number of mass condensates that give rise 
to caustic networks with higher complexity (Schneider & 
Weiss 1987; Blandford & Kovner 1988). Nevertheless, the 
models satisfactorily account for the symmetry and profile of 
the SAV features and serve as an “existence proof’ for the 
gravitational lensing hypothesis. 


5.2. Constraints on Lens Properties 


Having seen that simple lens models can produce the 
observed SAV features in J1415+1320, we now place 
constraints on the physical properties of the putative lens. 
Although several lens configurations are viable, fold caustics 
have generic properties that enable us to constrain the lens 
parameters almost independently of the precise lens 
configuration. 

We constrain the lens mass, M, by placing bounds on the 
Einstein radius, Og, that sets the natural length and angular scale 
for lensing (Blandford & Narayan 1992): 


0g = (4GM/c?D) ? 
e: 90(M /10? Mc)? (D/1 Gpe)! pas. (1) 


Here, G is the gravitational constant, and D = DiD,/Dy,, 
where D; and D, are the observer-lens and observer-source 
distance and Dj, is the lens-source distance. All distances are 
angular-diameter distances. Lensing causes multiple imaging 
with separations of up to zz20g. There is no evidence for 
multiple imaging in 15 GHz MOJAVE VLBI observations 
(Lister et al. 2009, 2016) at 23 epochs between 1995 and 2011. 
The observations had a resolution of about 1mas. By 
conservatively assuming that image splitting at 0.5 mas would 
be detectable, we place an upper bound on the Einstein radius: 
Og < 250 pas. Jets in AG usually exhibit moving features that 


propagate away from the core. If the jet is directed close to our 
line of sight the apparent speed can be superluminal (Begelman 
et al. 1984). In AGs at cosmological redshifts this probes lenses 
with Einstein radii ~100(v, /c)(t/1 year) uas. 

A lower bound on 0g can be placed using constraints on the 
sizes of synchrotron sources. For lensing to yield considerable 
magnification, the source must be significantly smaller than 0g. 
Fits to the precise profile of the observed SAV require 


p = 0,/0g < 0.03. (2) 


However, radio sources of a given flux density cannot have a 
size smaller than a threshold below which inverse Compton 
losses rapidly “cool” the relativistic plasma. This sets an upper 
limit of ~10'? K on the intrinsic brightness of an incoherent 
synchrotron source, though brightness temperatures of 
T, = 10? K have been observed in compact components in 
AG jets (Johnson et al. 2016), likely caused by relativistic 
Doppler boosting. If the source has a Gaussian brightness 
profile, the brightness temperature can be related to the 
observed source flux density by the Rayleigh— Jeans law: 


2 
DA nid. (3) 


where k is Boltzmann's constant, À is the wavelength, and ju is 
the lensing magnification. Since the Fermat potential near fold 
caustics has a generic functional form, the peak magnification 
near a fold also has a generic form (Schneider & Weiss 1992): 


lis m d opc. (4) 


So= pH 


where the constant of proportionality has been determined 
numerically (see Appendix B). Equations (1)-(4) can then be 
used to put a lower limit on the Einstein radius: 


So X 


ULT _. 
EZ 2.pl24skT, 


(5) 
The mass of the lens can then be determined from 0g for any 
given lensing geometry. 

In the case of J1415+1320, the peak-to-trough flux-density 
change for the U-shaped events is S, zz 0.5 Jy at 15 GHz and 
we assume p — 0.03. Figure 6 shows the lower bounds of the 
lens mass and Einstein radius for geometries where the radio 
source is located in the spiral and where it is a background 
radio source. Since the intrinsic sizes of the lensed source 
components are not precisely known, we plot the constraints 
corresponding to values of T, in the range of 107—10!^ K. 


5.3. The Location and Mass of the Lens 


Due to the ambiguity in the location of the radio source (see 
Section 4), we discuss both possibilities: (1) the radio source is 
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Figure 5. Multifrequency model light curves overplotted on data at 15, 37, and 234 GHz. The lens model and source trajectories are the same as that in the topmost 


row of Figure 4. 


located in the z = 0.247 spiral, and (ii) the radio source is an 
unrelated background object. In the former configuration, we 
can rule out all models in which the lens lies in our Galaxy; 
Galactic lenses consistent with the light curves of J1415+1320 
have a mass range consistent with optical microlensing 
surveys. Such surveys (Alcock et al. 2000) find lensing 
probabilities of about 107”, and crossing times of several days, 


1500 2000 2500 3000 


whereas we have found a crossing timescale of ~1 year. If both 
the source and the lens are in the spiral galaxy, the lens would 
need to be a 2106 Mo black hole. Since we do not find such 
objects in the spiral arms in our own Galaxy, we discount this 
possibility. The remaining possibility is that the lens lies 
between our Galaxy and the spiral host of the radio source. 
Here the requirements on the lens mass and projected density 


THE ASTROPHYSICAL JOURNAL, 845:89 (16pp), 2017 August 20 


Radio source in spiral 


Ty — 10? K, 0g — 194 uas 
Tp = 1013 K, Og = 61 uas 
T, = 10" K, dg = 19 pas 


Lens mass [Mo] 


10° 


10° 


104 


Projected density [Mo pc ?] 


10? 
1 2 3 4 5 6 7 


Di [105 pc], D, = 822 Mpc 


Vedantham et al. 


Background radio source 


10° 


10° 


104 


10° 


104 


1 2 3 4 5 6 
Di, [108 pc], D; = 822 Mpc 


Figure 6. Constraints on the minimum lens mass (top panels), and projected mass density (bottom panels). The left and right columns correspond to geometries where 
the radio source is in the putative spiral host, and where the radio source is an unrelated background object, respectively. The curves are plotted as a function of 


observer-lens distance, D;, and the source-lens distance, D;,, for the two geometries. 


are 210? Mo and Z10^ Mo pc~?. As discussed in Section 4.1, 
we obtained a deep optical spectrum of J1415+1320 with the 
Keck telescope and we did not find evidence of an intervening 
galaxy, although a faint dwarf galaxy could have evaded 
detection. In the scenario where the radio source is located in 
the spiral galaxy at z — 0.247 there are thus three awkward facts 
that must be accepted: (1) J14154-1320 is unique in being a BL 
Lac object and a CSO located in a spiral rather than in an 
elliptical galaxy, (ii) there is a 50 pc offset between the radio 
source and the optical nucleus, and (iii) the lens is an invisible 
intergalactic object, possibly a dwarf elliptical galaxy. 

In the alternative hypothesis the radio source is a background 
object and the foreground spiral galaxy provides a natural host 
for the lens. The absence of multiple imaging on arcsecond 
scales due to lensing by the spiral's bulge limits the redshift of 
the radio source to be z < 0.5 (Lamer et al. 1999). In addition, 
the angular size of interstellar scattering (ISS) in the spiral must 
be significantly smaller than the Einstein radius for ISS not to 
"quench" the lensing signatures. Non-detection of ISS at a 
wavelength of 18 cm (Perlman et al. 1996) shows that any 
scattering of the J1415+1320 emission by the spiral is low 
compared to several nuclear sightlines through our own Galaxy 
at low latitudes. However, ISS in our Galaxy is patchy, and 


little is known about ISS in other galaxies. In the background 
source scenario, there are two awkward facts that must be 
accepted: (1) the low probability of such an alignment, and (ii) 
the anomalously low ISS, but (1) is an a posteriori argument and 
Gi) is based on assumptions regarding ISS that cannot be 
substantiated by current observations. For these reasons, the 
background source argument is preferred, though not compel- 
ling. In any case, we see from Figure 6 that the similar 
constraints on the lens mass (>10° Mo) and projected density 
(210^ Mo pc?) apply to this geometry too. 


6. Discussion 


In this section, we consider the implications of the results 
obtained thus far and make falsifiable predictions to further test 
the hypothesis that gravitational lensing is the cause of SAV. 


6.1. Cosmological Density of Intermediate-mass Lenses 


The incidence of lensing may be used to constrain the 
cosmological fraction of matter in the putative lenses (Press & 
Gunn 1973). We adopt the Planck cosmological parameters: 
Ho —67.8kms ! Mpc !, Nm = 0.308, Q4 = 1 — On, Qk = 
0. Assuming a flat universe, this yields a critical density of 
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Figure 7. Fraction of lines of sight intercepted by the putative lens population 
(lensing optical depth) as a function of source redshift. The curves correspond 
to cases where different fractions of the cosmological matter density are in the 
lenses. The horizontal yellow line corresponds to the observed probability of 
intercepting a lens (based on just one sample, J1415+1320), and the shaded 
region is the corresponding 95% confidence interval assuming Poisson 
statistics. 


p, = 3H; / (81G) = 8.633 x 107? kg m ?. Let a fraction f 
of the matter be in the putative lenses of mass M. The 
comoving number density of lenses is then n = fp. /M. 
Consider a comoving volume element at redshift z; 
of dV (z) = (c/Ho)(1 + z} D? /E(z)dQ€dz, where E(z) = 
al (1 + zPQ,, + O4. There are, on average, ndV lenses in this 
volume element. Each lens presents a strong-lensing cross- 
section of n0?. The fractional sky area occupied by all lenses 
out to the source redshift is therefore 


1 
TG) = — [nav cone. (6) 
An 
Substituting for n, dV, and 0,, we get 
Z=Zs Q D 
T (Zs) = m f dz iu m ls 
Ho z=0 Cc DD, 
(1 + z? (7) 


x 
VC + zPO, + O4 


Figure 7 shows the expected lensing optical depth for 
various values of matter densities in the putative lenses 
Qı = fQm, as a function of source redshift. To estimate y, 
these curves must be compared to the observed probability of 
lensing—one lensed source among a sample of 981 objects, 
assuming J1415+1320 is unique in the sample. Because lensed 
objects have larger apparent brightness, they are overrepre- 
sented in any flux density-limited samples such as the CGRaBS 
sample. This magnification bias factor largely depends on the 
faint-end luminosity distribution of the specific source popula- 
tion under consideration, and the properties of the lens 
population. Both of these are difficult to establish, especially 
given our postulate of a new population of intermediate-mass 
lenses. We adopt a bias factor of 10 to derive an order-of- 
magnitude estimate. Accounting for this magnification bias, we 
find that one in 981 x 10 sightlines must be intercepted by the 
putative lenses. The errors on this lensing rate are large with 
just one sample, but they can be quantified by assuming 
Poisson statistics for the incidence of lensing. In Figure 7 we 
show the 95% confidence bounds on the lensing rate parameter 
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as a yellow box. Figure 7 shows that for a source redshift of 
z,— 0.5, a 95% bound of 10-^Q,, to 3 x 10-70, may be 
placed on the comoving density of the putative intermediate- 
mass lenses. Hence, if the lensing hypothesis is shown to be 
true, then a non-trivial fraction of the cosmological matter 
density must exist in intermediate-mass objects with projected 
densities of >10*M. pe ?. 


6.2. Falsifiable Predictions 


Our hypothesis that links SAVs to gravitational milli-lensing 
leads to several key predictions: 


(i) Intensive high-cadence long-term monitoring of J1415 
+1320 at multiple frequencies should detect further SAV 
events as new compact source components traverse 
the lens. 

(ii) The incidence of SAVs should increase with redshift (see 
Figure 8). 

(iii) During future SAV events, global millimeter-wave VLBI 
observations should be able to detect creation and 
annihilation of lensing images on ~1 year timescales 
(see also Gopal-Krishna & Subramanian 1996). The 
lensed images form at a typical separation of 
202 Z 100 uas which can be resolved by existing 
VLBI facilities. Bright image pairs will be closer than 
20g, however they will brighten and merge, or fade and 
separate in a manner that will be accessible to multi- 
epoch millimeter VLBI. 

(iv) The polarization position angle of the multiply lensed 
images should not change around caustic crossings. 

(v) Finally, if the recently discovered fast radio bursts 
(Lorimer et al. 2007) are at cosmological distances, then 
they must also be milli-lensed at a rate comparable to AG 
(Muñoz et al. 2016). In such cases, C10 ms time delays 
between lensed images of the burst might be detectable 
with upcoming surveys. 


6.3. SAV in the BL Lac Object AO 0235+164 


The prominence of the two SAV events in J1415+1320 in 
2009 and 2014 suggests the possibility that such signatures of 
radio milli-lensing may have been seen in other objects. For 
instance SAV may have been seen before, in the BL Lac object 
AO 02354-164. In the University of Michigan Radio Astron- 
omy Observatory low-frequency light curves of Aller et al. 
(2003, Figure 1), as well as the 22-90 GHz light curves of 
Nieppola et al. (2009, Figure 2), there is a symmetric crater 
U-shaped feature between 1994 and 1999, and the variability 
was achromatic from radio to optical frequencies (Webb et al. 
2000), although no optical observations are available for the 
first half of the crater feature. 

In Figure 8 we show the Metsáhovi 22 and 37 GHz light 
curves and the R and B band light curves (Webb et al. 2000) at 
high time resolution. Here the overall symmetry between 1994 
and 1999 is again clear, although there are also some departures 
from symmetry, such as might be caused by intrinsic variability 
in the radio source flux density and speed, and by the presence 
of multiple lensed components. The highly correlated varia- 
tions from radio to optical frequencies over the whole time 
interval could be produced by one radio component or several, 
moving behind a caustic network. Gravitational lensing has 
been invoked to explain various optical emission features of 
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Figure 8. The variability in AO 0235+164 at 22 and 37 GHz and B and R magnitudes (Webb et al. 2000; Nieppola et al. 2009). The roughly 5-year-long putative 


SAV is highlighted. 


AO 0235+164 (Stickel et al. 1988; Webb et al. 2000). The 
object has a redshift z=0.940 with intervening absorption 
systems at z=0.852 and 0.524 (Cohen et al. 1987), each of 
which could have one or several lenses associated with it. 
Detection of image splitting with high-resolution (<100 pas) 
global millimeter-wave VLBI observations is one avenue to 
ascertaining the presence or absence of milli-lensing in such 
SAV candidate sources. 

Itis worth emphasizing here that the two objects in which we 
find the clearest evidence of SAV, J14154-1320 and AO 0235 
+164, have both long been suspected of being subject to 
gravitational lensing (Stickel et al. 1988; Carilli et al. 1992; 
Stocke et al. 1992; Webb et al. 2000), and that our monitoring 
observations have led us, via the completely new and 
previously unsuspected phenomenon of SAV, to this same 
hypothesis. 


6.4. Outlook 


We have presented a rare new form of variability in AG that 
is time-symmetric and achromatic in nature. We propose that 
such variations are a result of gravitational lensing, wherein 
compact source components cross the magnification pattern of 
a gravitational lens at relativistic speeds. SAV may therefore 
provide a powerful new tool for investigating not only the cores 
of AG but also the putative lenses. These lenses have masses in 
the range 107—109 M..—a range that embraces intermediate- 
mass black holes, cores of globular clusters, dense molecular 
cloud cores, and compact dark matter halos (Carr 1994). The 
incidence of SAV in our sample suggests that a non-trivial 
fraction of the cosmological mass may be in such lenses, 
though meaningful constraints cannot be placed with just one 
or two examples. 

If SAV is confirmed to be due to gravitational milli-lensing it 
resolves the longstanding puzzles (Stickel et al. 1988; Carilli 
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et al. 1992; Stocke et al. 1992) over the high brightness and 
strong variability of the cores of J1415+1320 and AO 0235 
4-164, as well as the absence of a narrow emission-line region 
and the low thermal IR flux of J1415+1320. It would confirm 
the hypothesis posed some decades ago that a subset of 
variability in AG is due to gravitational lensing (Subramanian 
et al. 1985; Schneider & Weiss 1987; Gopal-Krishna & 
Subramanian 1991). In addition, it will bolster theories that 
posit that some members of the parent BL Lac population 
are inherently fainter and less variable (Ostriker & 
Vietri 1985, 1990), although these papers were suggesting 
permanent magnification and not what we are seeing in SAVs. 

In the 2014 SAV in J1415+1320 it appears that the 
millimeter radiation leads the centimeter radiation (Figure 1). 
This is not seen in the 2009 SAV, but there the millimeter 
sampling is sparse. Similarly, in AO 0235+164 there are some 
hints that the optical variations may precede the radio 
variations (see Figure 8). If our gravitational lensing hypothesis 
for SAV behavior is confirmed, and furthermore if these 
objects are not anomalous, they could have profound implica- 
tions for our interpretation of the relationship between the high- 
frequency and the low-frequency emission regions observed in 
AG jets. This behavior is readily interpretable if the moving 
features are outwardly propagating particle acceleration/ 
magnetic amplification fronts—shocks if the jet is plasma- 
dominated and nonlinear hydromagnetic waves if electro- 
magnetic field is dominant. The optically emitting electrons 
would cool close to the front and the centroid of their emission 
would cross any caustic ahead of the radio emission. Early 
shock-in-jet models (Marscher & Gear 1985) showed that 
millimeter variations lead centimeter variations in AG and 
would seem to fit this model. The alternative interpretation of 
multifrequency observations of AG jets is that the higher- 
frequency emission originates further back along the jet toward 
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the core, rather than further out along the jet away from the 
core, since it is well-known from VLBI observations that this is 
often the case at centimeter wavelengths in AG (Readhead 
et al. 1978; Readhead 1980). Whether this is the case in all AG 
and over a much wider frequency range is by no means clear. 
The additional resolution that SAV may provide under the 
lensing hypothesis could well then be the only direct way of 
exploring this interesting dichotomy other than via submilli- 
meter VLBI from space. 

The failure to associate dark matter with an elementary 
particle has revitalized interest in the idea that there may be a 
large population of primordial black holes (Carr et al. 2016). 
Further monitoring of relativistically moving radio source 
components in large samples of radio sources, combined with 
millimeter-wave VLBI imaging at <100 as scales may 
provide a new and sensitive approach to constrain this putative 
population. 

We end by noting that it is the time-symmetry and 
achromaticity of SAV that motivated the gravitational lensing 
hypothesis. If this hypothesis is falsified, then it would imply 
that SAV is intrinsic to the source, and pose an intriguing 
challenge to our understanding of relativistic jets. Attempts to 
reproduce intrinsic SAV features may provide critical insights 
into the emission geometry and mechanism in such jets. 
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Appendix A 
False Alarm Probability 


Lensing yields time-symmetric and achromatic light curves 
but they need not have the characteristic U-shape seen in J1415 
+1320. Here we have searched for U-shaped profiles in the 
OVRO data not necessarily to find milli-lensing candidates, but 
as a means to conduct a blind test of the null hypothesis: the 
SAVs in J1415+1320 are caused by stochastic fluctuations, 
where the "arrival" of U-shaped events is a Poisson process. 
Because the precise profile of the U-shaped events is not 
known a priori, automated searching for such events in our 
sample is non-trivial, and previous automated attempts 
(searches for ESEs) have only attained limited success (see for 
e.g., Lazio et al. 2001). We have embarked on an effort to 
develop algorithms for automated detection of SAVs, but here 
we employ human inspection (independently done by two co- 
authors) of the OVRO light curve to get an “order-of- 
magnitude" estimate of the false alarm rate. Figure 9 shows 
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the U-shaped events thus identified in the OVRO database. We 
note that, it is the repeating nature of U-shaped events in J1415 
+1320 that yield an extremely low false alarm-rate as 
computed in Section 3. Multiple U-shaped events are also 
seen in J0310--3814, but this object also shows a significantly 
larger amount of non-SAV variability than J1415--1320, and 
more importantly, it is unknown if the events in J03104-3814 
are achromatic. We also note that if shown to be chromatic, the 
U-shaped events in Figure 9 will be candidate ESEs (Fiedler 
et al. 1987). 


Appendix B 
Lens Modeling 


B.1. Rayshooting Simulation 


We employ the usual dimensionless units where all angles 
are expressed in units of the Einstein radius, 0g. Let the source 
and lens plane angular coordinates (two-dimensional) be 
y = [y y2] and x = [x x2], respectively. The lens equation 
relates the image-plane and source plane coordinates of light 
rays traveling from the source to the observer. For the single 
point-mass case, the lens equation is Schneider et al. (1992, 
chap. 8.1.2): 


yx (8) 


Ix 
where the second term on the right side is the gradient of the 
Fermat potential and is the angle by which the light rays are 
deflected by the lens. Consider a binary mass comprising two 
point-like objects of mass 1 and q located, respectively, at 
0—[0, 0] and d=[d, 0]. The lens equation then needs to 
include the contributions of both masses. 

The influence of an external smooth mass distribution can be 
parameterized by the convergence « and shear ^. The shear in 
particular accounts for any ellipticity in the external mass 
distribution. 


The total deflection is then (Schneider et al. 1992, 
chap. 8.2.2): 
x x—d Kd vy 0 
= x, (9) 
KP” x —aP | 0 «c. 


where “-” denotes matrix multiplication. We have chosen the 
principal axes for the ellipticity of the external mass to lie along 
the coordinate axes, so the off-diagonal elements in the last 
term of Equation (9) are zero. A positive y enhances ray 
deflection along the x, axes as compared to the x axes; this 
corresponds to an external mass distribution that is stretched 
along the x2 axis. Our models therefore make the implicit 
assumption that the binary axis lies along one of the cardinal 
axes of the external mass distribution's ellipticity. Relaxing this 
assumption extends the parameter space considerably and is not 
necessary to find adequate models. 

We used Equation (9) to simulate the magnification patterns 
using the technique of inverse-ray-shooting. We follow a 
large number of rays from the observer to the source plane. 
We “shot” about 10° rays on a regular grid on x. Using 
Equation (9), we compute the source plane location y of each 
ray. The rays are then binned on a suitable source plane grid. 
Because lensing preserves surface brightness, the magnification 
is the ratio of the image area to the source area on the sky. The 
magnification for each source bin is thus the ratio between the 
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Figure 9. OVRO 15 GHz light curves of seven objects in which U-shaped features were detected. The 10 features are indicated with arrows. 


number of rays in the bin and the number of rays in the absence 
of the lens deflection. The lens models used here are simple 
enough for the computation to be done on a desktop within 
reasonable compute-time («1 minute). 


B.2. Rationale for Choice of Lens Configuration 


We endeavored to generate a SRFD-FRSD “volcano” type 
light curve from a lens model that has the following properties: 
(i) the lenses should resemble commonplace astrophysical 
objects (point-like masses, elliptical halos, etc.), (ii) the lens 
model should be economical—simple models are favored; and 
(iii) a rectilinear source trajectory should yield the required 
SRFD-FRSD symmetry. 

We started with the simplest lens of all—an isolated point- 
mass lens (axially symmetric gravitational field). Such a lens 
has a point-like singularity in its source plane magnification 
pattern, and therefore yields a FRFD light curve, which cannot 
adequately describe the data. Elliptical perturbations to a 
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point-like lens typically yield a FRSD-SRFD (crater-type) light 
curve that is also not viable. Next, we considered binary lenses, 
which can yield a rich phenomenology of light curves for 
various values of the mass ratio q < 1 and separation d. We 
computed a library of magnification patterns by varying q and d 
and identified two scenarios that yield favorable magnification 
patterns: (i) an asymmetric (q « 1) binary with d ~ 1, and 
(ii) a binary with comparable masses (q X 1) with wider 
separation d = 1. The corresponding magnification patterns are 
shown in Figure 10, along with example source trajectories that 
give the SRFD-FRSD volcano-type light curve. For brevity, 
we will refer to the two lens models as “cupid’s bow" and 
“barbed arrow,” respectively. 


B.3. External Mass Sheet 


Having identified two promising candidate lens configura- 
tions, we then perturbed the lens parameters to account for 
(i) the relative timing of the four caustic crossings that together 
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Figure 10. Two binary lens configurations (top and bottom rows) that yield the required SRFD—FRSD parity. Top left: magnification pattern for a binary lens with 
q — 0.25, and d — 1.25 (the "cupid's bow"). Top right: magnification pattern for a source trajectory given by the black arrow. Bottom left and right: similar to the top 
panels, but for a lens with q = 0.01, and d = 0.8 (the “barbed arrow"). The top model can yield the two U-shaped parameters with a single source component, but the 
viable set of source trajectories is restrictive. The bottom model needs two source components to achieve two U-shaped events but admits a much larger range of 


source trajectories. 


yield the two U-shaped features, and (ii) the flux-density levels 
throughout the 2008—2016 period. These two constraints can be 
satisfied by the addition of a smooth background gravitational 
potential. Such a background potential could reasonably 
originate in the distributed mass of the galaxy that hosts the 
lens (Chang & Refsdal 1979). 

If k + |y| > 0.5, then the background mass has sufficient 
lensing strength to independently create multiple images on the 
macro scale (typically, arcseconds). Since no arcsecond-level 
image splitting has been observed in J1415--1320, we imposed 
the constraint & + || < 0.5. We again simulated a library of 
magnification patterns with this constraint and found that a 
large range of (x, y, q, d) parameter values can “fit” the light 
curves. The introduction of external mass is therefore required 
but is not very constraining. Two examples of favorable lens 
models that are perturbations of the cupid's bow and barbed 
arrow lenses are presented in Figure 4. 


B.4. Model Light Curve Fits to Multifrequency Data 


We obtained a model light curve fit to the OVRO 15 GHz 
data as follows. The data (Figure 1) show a linear decreasing 
trend in flux density over decade timescales, which may be due 
to variation in the unlensed source components. We subtracted 
the linear function of time, Su, (f) = —1.397 x 10-^t + 7.922, 
from the data, where flux density Sj is in Jy, and time t is MJD 


in days. Fitting model light curves to the data with the slope 
removed then requires the following parameters to be 
estimated: (i) the unlensed flux densities and angular sizes of 
the lensed components, (ii) the location and orientation of 
the source trajectory over the magnification pattern, (iii) the 
Einstein crossing time, fg, of the source components, and 
(iv) the flux density of source components that are not being 
lensed, assumed to be time-invariant. For simplicity, we 
assume that (i) all source components move at the same 
transverse speed along a north-south trajectory though the 
magnification pattern, and (ii) the source components have a 
bivariate Gaussian brightness profile with parameter [07, 07]. 
We fit the cupid’s bow lens with 0; = 6%, and the barbed arrow 
lens with 07 > 0%. Since the source size must be finite, we 
truncated the source brightness profile at +26,. 


Appendix C 
Possible Evidence for More Extended Fold Caustic Patterns 


The basic unit of symmetry we are exploring is the year-long 
U-shaped feature seen in 2009 and again in 2014 (Figure 1). 
However, there is possibly another symmetry evident in these 
data, namely the small peak near MJD 544664-1519, which is 
very close to the midpoint between the two peaks that straddle 
the small peak near MJD 54466--1519 in Figure 11, and are at 
MJD 54466+742.5+3.9 and MJD 54466+2295.7+1.8. These 
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Figure 11. OVRO 15 GHz data for the one-year interval centered on the light curve peak near MJD 544664-1519. The position of the midpoint between the two 
adjacent peaks (not visible here) and the corresponding uncertainty are shown by the vertical gray bar. This agrees with the centroid of the Gaussian fit to the central 
peak (dashed line) to remarkable precision (0.18%; see the text) and is possibly evidence for an extended caustic network. 


locations were determined by fitting low-order polynomials 
separately to the rising and falling parts of the light curves. The 
midpoint between these two adjacent peaks is therefore at MJD 
54466--1519.13-2.1. We have used the data shown in Figure 11 
to determine the position of the central peak by fitting a Gaussian 
to the central ~1 years of data, and we find that the centroid of 
the small peak lies at MJD 54466+1519.4+1.8. The difference 
of just 0.32-2.8 days between these two numbers may be 
compared with the relatively large interval of 1553 days between 
the adjacent two peaks. The above loc uncertainty amounts to 
only 0.1896 of the interval between the two peaks and may 
therefore indicate that the two U-shaped features observed in 
2008 and 2014 are related by a larger symmetry and that there is 
in fact a larger caustic network than we have assumed in the 
main text. Our simple lens models do not offer an explanation 
for this small peak at MJD 544664-1519, but given the 
symmetries that proliferate in gravitational lensing we believe 
that it should not be ignored at this early stage of exploring the 
SAV phenomenon. 


Appendix D 
Chance Alignment Probability 


Let there be Nr foreground galaxies to M; background radio 
sources and assume both are randomly distributed. We want to 
compute the expected number of instances where a background 
source and a foreground galaxy are projected (by chance) to 
within an angle sff from one another. The fraction of sky 
within a radial distance 6, from any of the N; foreground 
galaxies is 702. N;/(4x). The expected number of chance 
alignments is therefore MeN; 024 /4. If, however, the population 
of N, sources is drawn from a flux density-limited sample (such 
as the OVRO 40 m sample), then those in close alignment with 
a foreground galaxy will be overrepresented due to magnifica- 
tion by the putative millilenses. We assume an average milli- 
lensing magnification of u = 10 (see Figure 6). For an intrinsic 
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source-count power-law with slope of 1.5 (Euclidean value), the 
magnification bias factor is B, = 10'° ~ 30. For the OVRO 
survey, we have N, = 981, and the required alignment for 
J1415+1320 is Op = 13 mas. Based on the fits to the Press 
—Schechter function of Tomczak et al. (2014), we estimate that 
there must be about Ny = 2 x 108 halos above a mass of 10? Mo 
up to redshift z — 0.5 (putative redshift of J14154-1320). Hence 
the expected number of radio sources in our sample with the 
required alignment is B,N; N, 02. / 4 zz 0.006. While this is a low 
value, we remind the reader that this is a case of a posteriori 
statistics, so the argument against the background source 
hypothesis is not compelling. If halos of 108 and 107 Mo are 
also considered as viable foreground systems, then the corresp- 
onding numbers of chance alignments are 0.03 and 0.2, 
respectively, although it is unclear if halos in this mass range 
will host the putative lenses. 
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